Home | | Latest | About | Random
# Week 1 Tuesday
## Snow fall discussion and a discrete simulation.
For our snow fall problem
>*Snow falls before noon,*
>*falling at a constant rate.*
>*Jack starts shoveling snow at noon,*
>*shoveling at a constant rate.*
>*The first hour Jack shovels one mile,*
>*but the next hour, only half a mile.*
>*When did snow fall begin?*
Note that from the first half
>*Snow falls before noon,*
>*falling at a constant rate.*
>*Jack starts shoveling snow at noon,*
>*shoveling at a constant rate.*
We model it as a differential equation for Jack's position $x(t)$ as a function of time $t$: $$
x'(t) = \frac{c}{w} \frac{1}{at+b}
$$where $a,b,w,c$ all some positive constants. This is done by imposing a lot of assumptions. Also in our problem, we modeled that the time the snow starts falling is $t^{\ast} = -b /a$.
But now that we have a differential equation, we can solve it (via FTC), and obtain: $$
x(t) = \frac{c}{w a} \ln (at+b) + k
$$
for some constant $k$.
Now with the condition of Jack starts at $x(0)=0$ at $t = 0$, we determined $k$ and get $$
x(t) = \frac{c}{w a} \ln\left( \frac{a}{b}t +1 \right)
$$after some algebra. Notice here $x(t)$ is a logarithmic function in $t$.
And finally with the second half of the problem:
>*The first hour Jack shovels one mile,*
>*but the next hour, only half a mile.*
which translates to $x(1) = 1$ and $x(2) = 3 / 2$, we are able to deduce the time $t^{\ast}$ when the snow starts, by solving for the ratio of $a / b$.
We might scrutinize our model, is it descriptive enough?
Let's try to see if we can corroborate our model by doing a simplified discrete simulation of this snow fall problem.
Let us say the road is chopped into discrete segments, and on each segment we record the amount of snow it has on it. Initially, the road has no snow on it. We can represent this as an array of numbers. In python,
```
road = [0] * 100 # represent a road of length 100 with initially no snow on it
```
Now, we simulate snowing on the road by adding some amount of snow per time step. We can define a function to do that
```
snow_rate = 1 # amount of snow to add on the road per time step
def snow(road): # adds snow onto the road
for i in range(len(road)):
road[i] += snow_rate
return road
```
Next, we simulate shoveling, at some rate
```
shovel_rate = 3 # amount of snow we remove in each time step
def shovel(road,position): # remove snow on the road at given position
road[position] -= shovel_rate
return road
```
Now let us record the shoveler's position at each time step in an array, and update the position to increment if the snow at the current position is no longer positive:
```
positions = [0] # initial position
for t in range(1000): # simulate 1000 time steps
current_position = positions[-1]
road = snow(road)
road = shovel(road, current_position)
if road[current_position] <= 0:
positions.append(current_position + 1)
else:
positions.append(current_position)
```
We can print out the shoveler's latest position after simulating some number of time steps:
```
print(positions[-1])
```
Also, since we have a list of positions at each time step, we can plot this out
```
import matplotlib.pyplot as plt
plt.plot(positions)
plt.show()
```
And we get a graph that looks logarithmic! This corroborates with our expectation of the behavior of $x(t)$ with our differential equation model.
Note: You can adjust various parameters above to see how the simulation changes. Also, we didn't account for the time the snow falls before Jack starts shoveling, but you can add that in as well if you like. The point is -- we often do computer simulation to check whether our modeling is reasonable.
## Separation of variables.
- Recall that the differential equation $y' = \frac{1}{1+x^{2}}$ for $y=y(x)$ can be solved via FTC, but we get stuck for $y' = \frac{1}{1+y^{2}}$ for $y=y(x)$. Since naively applying FTC gives an **integral equation** $$
y(x)= y(0) + \int_{t=0}^{x} \frac{1}{1+y(x)^{2}} dx .
$$While this equation is true, it doesn't tell us what $y(x)$ is, since its expressing depends on an integral that involves $y(x)$ itself.
- However, all is not lost. This equation is what we called **separable equation**. Let us write our function as $y=y(x)$, and its derivative as $y'=\frac{dy}{dx}$. A first-order differential equation is called separable if it has the form $$
\frac{dy}{dx} = A(x)B(y)
$$where $A(x)$ is an expression of $x$, and $B(y)$ is an expression of $y$. That is, $\frac{dy}{dx}$ is a product of an expression purely in $x$, and an expression purely in $y$.
- If we "pretend" $dy / dx$ is a fraction, we can separate the variables as follows: $$
\frac{dy}{B(x)} = A(x)dx
$$Since the context of $dx$ or $dy$ that we know of besides in Leibniz's derivative expression is an integral, we integrate both sides: $$
\int \frac{dy}{B(y)} = \int A(x) dx.
$$The point is, if we have a separable equation $y' = A(x)B(y)$, it must satisfy the integral relation $\int \frac{dy}{B(y)} = \int A(x)dx$.
- Proof. Starting with$\int \frac{dy}{B(y)}$, we can perform a substitution of $y=y(x)$, where $dy = y'(x)dx$. So we get $\int \frac{y'(x)dx}{B(y)}$. But we know $y'=A(x)B(y)$, so this gives $\int \frac{A(x)B(y)dx}{B(y)} = \int A(x)dx$ as claimed. $\blacksquare$
- Remark. With this, we know the manipulation to get to the integral expression is legitimate. However, often to get there, we do use the trick of "treating $dy / dx$" as a fraction. As it turns out, we can properly define what $dy$ or $dx$ is as a mathematical object called **differential forms**. This allows us to expression first-order ODEs in the form of $M(x,y)dx + N(x,y)dy = 0$.
- For example, $y' = \frac{1}{1+y^{2}}$ has differential form expression $dx - (1+y^{2})dy = 0$.
- So back to our differential equation for $y=y(x)$, $y' = \frac{1}{1+y^{2}}$, we have $$
\frac{dy}{dx} = \frac{1}{1+y^{2}} \implies (1+y^{2})dy = dx
$$so integrating both sides give $$
y+\frac{1}{3}y^{3} + C_{1} = x+C_{2}
$$for some constants $C_{1}, C_{2}$. Since they are arbitrary constants, we can combine them into one, and obtain an expression $$
y + \frac{1}{3}y^{3}+C = x.
$$
- Note that we obtained an **implicit expression** for $y$. That is, $y=y(x)$ is a function that satisfies above equation, but we didn't explicitly express $y$ as a function of $x$. This is called an **implicit solution** and often it is the next best thing if we don't have a nice explicit way of writing our solution. Not always we can find explicit solutions easily.
- Also note, we have an arbitrary constant $C$. This means for each $C$, we get a different solution. So we didn't just get one solution, bur actually a **1-parameter family of solutions**, where $C$ is the parameter. This is a collection of solutions, one per choice of $C$. Try plotting them in Desmos. They are just side-ways cubics.
- If we are further given some initial conditions, say $y(2) =3$, then we can substitute to determine which $C$ it is. In this case, we would have $3 + \frac{1}{3}3^{3} + C =2$, so $C=-10$. So an implicit solution to $y' = \frac{1}{1+y^{2}}$, $y(2) =3$ is $y+\frac{1}{3}y^{3}-10 = x$.
- Remark. Not every differential equation is separable. For example: $y' = x+y$ for $y=y(x)$. This we require another technique. (Do you recognize what kind of DE it is? That's right, it's linear!)
- Another example. Find an explicit solution to $y' = x+y+xy+1$ for $y=y(x)$ satisfying $y(2)=1$.
- Answer. Note we have $\frac{dy}{dx} = (1+x)(1+y)$, which is a separable equation. So we set up $$
\int \frac{1}{1+y}dy = \int(1+x)dx \implies \ln|1+y|= x+\frac{x^{2}}{2} + C
$$Note how I combined the two constants of integration as one $C$. Now we can actually explicitly solve for $y$ here. Exponentiating both sides gives $$
|1+y| = e^{x + x^{2}/2+C} = e^{C}e^{x+x^{2}/2}
$$Now we are going to perform something called **abuse of constant**: Since $e^{C}$ is just a constant, we can relabel it as a constant $C_{1}$. Now, as there is no other constants around, we might as well recycle the symbol $C$ for this. (If this bothers you, you can just call it $C_{1}$.) So we have $$
|1+y| \overset{!}= C e^{x + x^{2}/2}
$$(**I use exclamation mark ! to indicate we abused the constants**) Now to take away the absolute values sign, we know that we get $1+y = \pm Ce^{x+x^{2}/2}$. But wait, $\pm C$ is just another constant. So by abusing the constant again, we get $$
1+y \overset{!}= C e^{x+x^{2}/2}.
$$Finally we have $y = C e^{x+x^{2}/2}-1$. This is a 1-parameter family of solutions to $y'=x+y+xy+1$. They are expressed explicitly now. Finally, to find the one satisfying $y(2)=1$ we just substitute and solve for $C$: $1 = C e^{4}-1$, so $C=2 / e^{4}$. This gives the explicit solution $y = \frac{2}{e^{4}}e^{x+x^{2}/2} - 1$.
- By the way, we can also rewrite $y' = x+y+xy+1$ as $y'-(1+x)y=x+1$, which is a linear DE. Later when we see how to solve first order linear ODEs you can check if we get the same answers.
## M is for ...
On a cold night in Fargo, a dead body was found in the forest!
When the police arrived, they measured the body's temperature to be $50$ degrees at 10 pm.
And at 11 pm, the body is measured to be $40$ degrees.
The weather reports that the whole day Fargo had a temperature of $30$ degrees.
When was the person murdered?
Let us assume that a living person has body temperature of $98$ degrees. And as they die they stop maintaining this temperature.
Let us denote $T(t)$ to be the temperature in degrees of the body at time $t$. And let us take $t = 0$ to mean 10 pm, and $t$ measured in hours.
To model this problem we have **Newton's law of cooling and heating** (which, he also found empirically -- actually they are still trying to define what temperature means back then). But it states the following:
> The rate of change of an object's temperature is proportional to the difference between the object and the environment.
So if we denote $T(t)$ as the temperature of the object at time $t$, and $T_{env}$ as the environment's temperature, we have $$
\frac{dT}{dt} = k (T-T_{env})
$$for some proportionality constant $k$. For simplicity sake, $T_{env}$ is also a constant. Here we can do a sanity check of the sign of $k$. Written this way, if $T > T_{env}$, then physical intuition says $dT / dt < 0$, so $k$ would be negative; similarly if $T < T_{env}$, then $dT / dt > 0$, and still $k < 0$. Ok, let us continue.
Aha! This is a separable equation! So we have differential form $$
\frac{dT}{T-T_{env}} = k dt
$$which we can integrate to get $$
\ln|T-T_{env}| = kt + C
$$after combining constants to $C$. Exponentiating and abusing constants as before, we get $$
|T-T_{env}| = e^{kt + C} \overset{!}= Ce^{kt}
$$so $$
T-T_{env} \overset{!}= Ce^{kt}
$$This gives a 1-parameter family of curves, with parameter $C$.
We can apply all the information we have now: Since $T_{env} = 30$, so $$
T=30 + Ce^{kt}
$$
At $t = 0$, which is 10 pm, we have $T=50$. So $$
50 = 30 + Ce^{k (0)}=30 + C \implies C = 20
$$
At $t=1$, which is 11 pm, we have $T = 40$. So $$
40 = 30 + 20e^{k} \implies e^{k} = \frac{1}{2} \text{ or } k=\ln\left( \frac{1}{2} \right)
$$
Hence $$
T (t) = 30 + 20 e^{t\ln(1 / 2)} = 30 + 20 \left( \frac{1}{2} \right)^{t}.
$$
Now we can find when the person was killed, if living temperature is $98$: Solve for $t$ where $$
98 = 30 + 20 \left( \frac{1}{2} \right)^{t}
$$which gives $$
t = \frac{\ln\frac{68}{20}}{\ln \frac{1}{2}} \approx -1.7655 \text{(hrs)}
$$
So that is $\approx 10 -1.7655 = 8.2345$ pm, which, to the nearest minute, is about $8:14$ pm!
## First-order linear ODEs.
A first-order linear ODE has the form $$
a_{1}(x) y'(x) + a_{0}(x)y(x) = f(x)
$$
If $a_{1}(x)$ is not zero, then we can make it **monic**, namely the leading coefficient is 1. We will be concerned with **monic first-order linear ODEs** that has the form $$
y'(x) + P(x) y(x) = Q(x)
$$for unknown function $y(x)$, given $P(x)$ and $Q(x)$.
Looking at the expression we are inspired by product rule. Note if $\mu(x)$ and $y(x)$ are two functions of $x$, then the derivative of their product is $$
[\mu(x)y(x)]' = \mu'(x)y(x) + \mu(x) y'(x).
$$
So let us take a yet-to-be determined function $\mu(x)$, called an **integrating factor**, and we multiply our monic first-order ODE, and get $$
\mu(x)y'(x) + \mu(x)P(x) y(x) = \mu(x) Q(x).
$$
Now if we impose the condition that $\mu'(x) = \mu(x)P(x)$, then we can recognize the left-hand side above as the derivative of the product, and get $$
[\mu(x)y(x)]' = \mu(x)Q(x).
$$Then by FTC, we see that $$
\mu(x)y(x) = \int \mu(x)Q(x)dx + C
$$or, solving for $y(x)$, we get $$
y(x) = \frac{1}{\mu(x)} \left[ \int \mu(x)Q(x)dx + C \right]
$$This works provided that $\mu(x)$ satisfies $\mu'(x) = \mu(x)P(x)$. But wait, this is separable and we can solve for $\mu$! Note $$
\frac{d\mu}{dx} = \mu(x) P(x) \implies \frac{d\mu}{\mu} = P(x)dx.
$$So $\ln|\mu| = \int P(x)dx$. By abusing constants we get that $\mu(x) = C_{1}e^{\int P(x)dx}$, where $C_{1}$ some constant. But! If you look at the expression of $y(x)$ we obtained, this $C_{1}$ constant can be canceled as $\mu$ shows up in both the numerator and denominator; and $C_{1}$ can be combined with $C$ if we distribute. So it suffices to take $$
\mu(x) = e^{\int P(x) dx}.
$$
So a summary. For the monic first-order linear DE of the form $y'(x) + P(x)y(x) = Q(x)$:
- Write down the integrating factor $\mu(x) = e^{\int P(x) dx}$.
- Write down the solution $y(x) = \frac{1}{\mu(x)} \left[ \int \mu(x)Q(x)dx + C\right]$.
Let us see an example. Find a 1-parameter family of solutions to $y' = x + y$, for $y=y(x)$.
Let us rearrange it as $y'-y=x$. So we recognize $P(x)=-1$ and $Q(x) = x$. Our integrating factor is then $\mu(x) = e^{\int-1 dx} = e^{-x}$. So we can write down the solution of the form $$
\begin{align*}
y(x) & = \frac{1}{\mu(x)} \left[ \int\mu(x) Q(x) dx + C\right] \\
& = e^{x} \left[ \int e^{-x}x dx + C \right]
\end{align*}
$$We can solve for the integral there by integration by parts, where $$
\begin{align*}
\int e^{-x}x dx &= \int x d(-e^{-x}) \\
&= -xe^{-x} - \int -e^{-x} dx \\
&= -x e^{-x} -e^{-x} +C
\end{align*}
$$
So as our 1-parameter family of solutions, we get $$
y(x) = e^{x} [-x e^{-x} - e^{-x} + C]
$$ or $$
y(x) = -x-1+Ce^{x}.
$$
We can verify this by substituting this back and check if $y' = x+y$. Note that $$
y' = -1 + Ce^{x}
$$and $$
x+y = -1+C e^{x},
$$which they agree!